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I.  INTRODUCTION 

A.  BACKGROUND 

The  Strategic  Defense  Initiative  Organization  (SDIO)  has  developed  a  goal 
to  vigorously  research  and  develop  technology  that  could  help  to  eliminate  the 
threat  of  ballistic  missiles  and  provide  increased  U.S.  and  allied  security.  By 
deploying  a  three-part  phased  ballistic  missile  defense  system,  incremental 
strategic  benefits  can  be  realized  while  preparing  the  way  for  the  next  phase. 

The  first  phase  would  reduce  confidence  of  planners  initiating  a  nuclear 
attack  against  the  U.S.  by  not  allowing  them  to  predict  the  outcome  of  a 
ballistic  missile  attack.  The  second  phase  would  negate  the  potential  threat 
government’s  or  hostile  organization’s  ability  to  destroy  many  of  the  U.S. 
strategic  targets,  and  the  third  phase  would  eliminate  the  threat  posed  by 
ballistic  missiles  entirely  (Udall,  1988). 

The  first  phase  system  being  proposed  by  SDIO  includes  ground  and  space 
based  BMD  consisting  of: 

•  Spaced-based  hit-to-kill  vehicles  for  attacking  missile  boosters  and 
post-boost  vehicles. 

•  Ground  based  rockets  designed  to  intercept  warheads  as  they  reenter 
the  atmosphere. 

These  spaced  based  hit-to-kill  interceptors  (SBI’s)  would  be  arranged  in 
a  constellation  of  several  hundred  satellites  at  several  hundred  kilometers 
altitude  above  the  earth.  A  constellation  of  satellites  is  an  organized  collection 
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of  satellites  in  related  orbits.  Each  satellite  would  have  the  capability  to  detect 
ballistic  or  anti-satellites  missile  launches  by  observing  the  hot  rocket  plumes. 

Once  a  ballistic  missile  has  been  detected,  the  SBI  would  be  able  to  track 
the  booster  and  pass  this  information  to  the  rest  of  the  constellation.  This 
information  sharing  has  two  uses: 

•  Due  to  positional  constraints,  the  SBI  that  is  tracking  a  booster  may 
not  be  the  SBI  that  has  the  best  shot  at  killing  the  booster.  Tracking 
information  can  be  passed  to  the  SBI  with  the  best  shot 
(Compaterto,  1991). 

•  Line  of  sight  laser  communication  will  be  used  to  minimize  jamming. 
This  requires  sequential  message  transmission  from  one  SBI  to  another 
SBI  and  could  result  in  large  queues  of  messages  being  formed 
(Comparetto,  1991). 

The  use  of  a  constellation  of  orbiting  SBIs  to  identify,  track  and  engage 
thrusting  bodies  has  numerous  advantages.  The  system  can  be  made  to 
operate  autonomously,  provide  world  wide  coverage,  and  it  is  flexible  to 
changing  political  situations. 

A  constellation  of  SBIs  will  consist  of  hundreds  of  platforms,  each  with 
identical  capability.  The  large  number  of  platforms  will  result  in  multiple 
coverage  of  any  given  area.  Given  a  booster  launch,  more  than  one  SBI  will 
observe  the  launch  and  commence  tracking  that  booster.  Due  to  different 
observation  angles,  the  position  of  the  sun,  and  the  individual  sensor  systems 
*  themselves,  the  tracking  quality  of  these  SBIs  will  be  variable.  Some  tracking 
information  will  be  better  than  others.  Information  of  high  quality  should  be 
communicated,  while  poor  quality  information  should  not  be  passed  on  to  other 
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SBIs.  This  type  of  pruning  will  reduce  time  in  queues,  decrease  the  time  for 
information  being  transmitted  throughout  the  constellation,  and  will  allow  the 
SBI  with  the  best  shot  a  higher  probability  of  kill  (Comparetto,  1991). 

The  software  developed  to  simulate  a  constellation  of  SBIs  is  one  of  the 
Strategic  Defense  System  (SDS)  Simulators.  The  template  based  tracking 
algoritbm  is  a  function  of  this  simulator  and  tracks  booster  and  ballistic  bodies 
using  a  single  optical  tracker.  It  is  this  papers  goal  to: 

•  Describe  the  functionalities  of  the  template 
based  tracking  algorithm  and  how  it  works. 

•  Determine  if  there  is  a  reliable  way  to  measure  the 
quality  of  a  track  message  generated  by  the  tracking 
algorithm. 

•  Make  appropriate  conclusions  and  recommendation  to 
improve  the  template  based  tracking  algorithm. 

Chapter  II  will  give  a  brief  description  of  the  system  simulator  and  how 
the  track  algorithm  works.  Chapter  III  will  describe  a  measure  of  quality  for 
a  track  message  and  test  this  measure  using  the  tracking  algorithm. 

Chapter  IV  will  show  the  effect  of  changing  the  azimuth  and  elevation  variance 
of  the  sensor  c  *1  the  tracking  algorithms  error.  Chapter  V  will  give  the 
conclusions  and  recommendations. 
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H.  DESCRIPTION  OF  THE  SDS  SIMULATOR 


A.  OVERVIEW 

The  system  simulator  uses  the  template  based  tracking  algorithm  to 
estimate  the  launch  parameter  state  vector  and  the  launch  parameter  variance- 
covariance  matrix  of  ballistic  targets  based  on  data  from  a  single  optical  sensor. 
This  problem  generally  has  no  unique  solution  unless  additional  information  is 
used.  Multiple  trajectories  can  be  constructed  through  the  same  measured 
angle  data.  Constraints  are  required  to  make  the  problem  solvable.  The 
approach  taken  by  the  tracking  algorithm  is  to  utilize  a  priori  trajectories  in 
the  form  of  downrange  and  altitude  templates  which  are  specific  for  a  given 
booster  type.  These  templates  consist  of  a  family  of  curves,  each  one 
representative  of  a  flight  profile,  and  in  the  ensemble  encompass  the  flight 
envelope  of  a  given  booster  type  (Figure  2-1,  Rasmussen,  1989). 


Figure  2-1.  Multiple  trajectories  constructed  from  measured  angle  data. 
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B.  COORDINATE  SYSTEM 


The  template  based  tracking  algorithm  uses  five  different  coordinate 
systems.  The  first  is  the  Earth  Centered  Inertial  (ECI)  coordinate  system. 
The  x  axis  is  contained  in  the  earth’s  equatorial  plan  and  is  directed  through 
the  prime  meridian,  the  z  axis  is  directed  through  the  north  pole,  and  the  y 
axis  is  perpendicular  to  both  completing  the  right-handed  system.  The  second 
system  is  the  Earth  Fixed  (EF)  coordinate  system  which  is  aligned  to  the  ECI 
system  at  epoch,  but  whose  x  and  y  axes  rotate  at  the  rate  of  the  earth.  The 
third  coordinate  system  is  the  geographic  coordinates,  consisting  of  latitude, 
longitude,  altitude  and  launch  azimuth.  The  fourth  coordinate  system  is  the 
local  launch  coordinate  system,  where  the  axis  is  contained  in  the  local  tangent 
plane  and  is  directed  along  the  launch  azimuth  of  the  target,  the  z  axis  is 
directed  towards  local  zenith,  and  the  y  axis  completes  the  right  handed 
system.  The  fifth  coordinate  system  is  the  sensor  local  coordinate  system 
where  the  x  axis  is  contained  in  the  local  tangent  plane  and  is  directed  locally 
north,  the  y  axis  is  contained  in  the  local  tangent  plane  and  directed  locally 
east,  and  the  z  axis  is  directed  to  local  zenith  (Rasmussen,  1989).  This  is  a  left 
handed  system  (Figure  2-2,  Rasmussen,  1989). 

C.  FUNCTIONALITY 

The  problem  of  estimating  the  trajectory  of  a  thrusting  target  using  only 
the  angle  of  measurement  from  a  single  optical  sensor  is  not  well  defined. 
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Multiple  trajectories  can  be  constructed  through  the  measured  angle  data.  To 
solve  this  problem,  trajectory  information  in  the  form  of  downrange  and 
altitude  templates,  which  are  specified  per  booster  type,  are  needed 

The  trajectory  templates  are  given  as  the  downrange  and  altitude  of  a 
booster  as  a  function  of  time  for  that  booster  type.  Various  pitch  profiles  are 
included  to  take  into  account  lofted  or  depressed  trajectories.  These  curves  are 
used  as  constraints  by  the  template  based  tracking  algorithm.  It  is  an  implicit 
assumption  that  any  particular  booster  trajectory  may  be  approximated  as  a 
linear  combination  of  these  a  priori  trajectory  templates.  Since  the  trajectory 
templates  encompass  the  total  dynamic  behavior  of  a  booster’s  trajectory,  any 
large  deviation  from  the  nominal  shape  of  the  altitude  and  downrange 
templates  serves  to  degrade  the  algorithm  performance  (Rasmussen,  1989). 

As  an  internal  function  of  the  tracking  algorithm,  the  trajectory  templates 
are  represented  by  bicubic  splines.  Cubic  splines  are  constructed  to  fit  the 
altitude  and  downrange  templates  for  each  flight  profile  independently.  At  any 
given  iteration  in  the  launch  parameter  estimation,  the  altitude  and  downrange 
are  interpolated  from  the  cubic  splines  for  each  flight  profile  based  upon  the 
estimated  time  of  flight.  A  cubic  spline  is  then  constructed  at  the  estimated 
time  of  flight  across  the  flight  profiles,  and  an  interpolated  value  is  obtained 
based  upon  the  current  estimate  of  the  pitch  parameter.  In  order  to  estimate 
the  partial  derivatives  necessary  to  determine  the  gradient,  altitude  rate, 
altitude  acceleration,  downrange  rate,  and  the  downrange  acceleration  are 
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evaluated  by  differentiating  the  cubic  spline  polynomials  for  each  flight  profile. 
Cubic  splines  are  then  constructed  across  the  rate  and  acceleration  points  and 
evaluated  at  the  current  pitch  parameter  estimate. 

D.  LAUNCH  PARAMETER  ESTIMATION 

The  template  based  tracking  algorithm  estimates  the  launch  parameter 
using  an  iterative  batch  least-square  algorithm.  If  X  is  a  six  dimension  state 
vector  of  the  launch  parameters  and  Z  is  the  azimuth  and  elevation 
measurement,  then  the  relationship  exists: 

Z  =  h(X)  +  v 

where  the  h(X)  is  a  function  dependent  on  the  a  priori  booster  information  and 
v  is  multivariate  normal  error  with  mean  zero  and  a  given  variance-covariance 
R  .  Taking  a  first  order  Taylor  expansion  about  X0,  the  fixed  state  vector  of 
the  best  known  launch  parameter  gives: 

Z  -  h(XJ  +  (X-XJH(XJ  +  v 
where  H  is  the  derivative  of  h  evaluated  at  X,.  If  we  let: 

fiZ  =  Z-h(XQ),  5X  =  (X-XJ 

then  it  follows  that: 

fiZ  »  H(XJfiX  +  v 

From  this  equation  using  Least  Squares,  fiX,  which  is  the  difference  between 
the  best  known  launch  parameter  state  vector  and  the  new  estimate  for  the 
state  vector,  can  be  estimated. 


Let  the  variance-covariance  matrix  of  v  be  defined  as  R  From  least 


squares  (Mendenhall,  1989)  an  estimator  of  5X  is: 

fix’  =  ((R1/*H)T(R1/*H) ■1/I(R*1/1H)T(R-1/1fiZ) 
=  (HTR'1H)1HTR'1fiZ 

where  R  =  R1/5R1/J  and  H  =  H(XJ.  This  only  works  if: 

E[5Z]  =  H(XJ5X 


but  in  fact: 


E[5Z]  =  h(X)  -  h(XJ  -  H(X0)5X 
An  estimator  of  X  is  therefore: 


X,  =  X„  +  fix’ 

Therefore  X,  is  the  estimate  of  the  launch  parameters.  Because  a  linear 
approximation  of  h  about  X„  is  used,  replacing  X,  for  X„  and  solving  the  least 
squares  problem  iteratively,  the  estimator  of  the  launch  parameter  at  the  nth 
iteration  is: 


x„  =  x,.,  +  fixv, 

The  position  of  the  booster  in  Earth  Centered  Inertial  (ECI)  is  a  function  of  the 
launch  parameters,  based  on  time,  downrange  and  altitude  given  by: 

X^j  =  Z(-w(t-t,1>))(Xo(0,X,ho)  +  TC^.X.oJX^Ch^r,^)) 


where: 

Z(-co(t-t^,))  is  the  transformation  matrix  from  earth 
fixed  coordinates  to  ECI  coordinates 

XD(^,X,h0)  is  the  launch  position  vector  expressed  in 
earth  fixed  coordinates 


9 


T(4>,\,a)  is  the  transformation  matrix  from  local 
launch  coordinates  to  earth  fixed  coordinates 

X~.i(h„.T,fl)  is  the  current  booster  position 

alone  the  trajectory  template  expressed  in 
local  launch  coordinates 


The  covariance  matrix  associated  with  the  launch  parameter  estimate  can  be 
approximated  in  ECI  coordinates  by  the  following  linearized  equation: 

P  *  (5X/6Y)C(5X/6Y)t 


where: 

P  is  the  covariance  matrix  of  the  ECI  state  vector 

C  is  the  covariance  matrix  of  the  launch  parameter 
vector  estimate 

X  is  the  ECI  state  vector 

Y  is  the  launch  parameter  vector 

5X/5Y  is  the  Jacobian  of  the  transformation  form 
launch  parameters  to  ECI  coordinates 

The  template  based  tracking  algorithm  computes  a  variety  of  output  including: 
estimated  ECI  position,  velocity,  acceleration,  covariance  matrix,  launch 
parameters  and  launch  parameter  covariance  matrix  for  the  booster. 
Additionally  the  true  ECI  position  of  the  booster,  ECI  position  of  the  sensor 
platform,  and  other  parameters  are  available  (Figure  2-3,  Rasmussen  1989). 
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Figure  2-3  Estimated  errors  in  the  tracking  algorithm. 
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m.  MEASURE  OF  ERROR 


A.  COVARIANCE  MATRIX 

It  was  believed  that  the  best  measure  of  error  in  the  system  would  be  the 
estimated  ECI  variance-covariance  matrix  itself.  Initially,  it  was  hypothesized 
that  large  differences  between  the  estimated  position  and  the  actual  position 
of  the  thrusting  body  would  be  reflected  in  the  estimated  launch  parameter 
covariance  matrix  and  the  estimated  ECI  variance-covariance  matrix  .  The 
launch  parameter  is  expressed  in  geographical  coordinates,  time,  altitude,  pitch 
and  azimuth.  The  ECI  covariance  matrix  is  expressed  in  kilometers,  a  natural 
measure  of  distance  and  error. 

The  version  of  the  tracking  algorithm  used  was  the  mini  testbed,  which 
was  developed  as  an  analysis  tool.  As  such,  it  is  flexible  in  options  such  as  scan 
rate  (number  of  samples  taken  of  the  simulation  run),  the  initial  position  of 
both  the  booster  and  the  sensor  platform,  the  booster  type,  booster  data,  and 
simulation  run  time.  All  runs  used  simulation  runs  of  100  seconds,  and  the 
same  initial  positions  for  the  sensor  and  booster.  A  scan  rate  of  2  seconds  was 
used. 

The  tracking  algorithm  relies  on  the  assumption  that  the  errors  are 
multivariate  normal.  Initial  runs  of  the  tracking  algorithms  were  conducted  to 
discern  the  distribution  of  the  error.  Plots  in  the  ECI  X,  Y,  and  Z  coordinate 


plans  were  made  by  taking  the  difference  in  the  estimated  X  (  or  Y  or  Z) 
position  and  the  real  X  (or  Y  or  Z)  position.  The  shape  of  the  distributions  was 
bell  shaped  and  generally  centered  around  0.  The  distribution  is  not  spherical. 
(Figure  3-4). 

The  assumption  that  the  distribution  is  multivariate  normal  with  mean  0  and 
an  unknown  variance,  is  not  unreasonable  (Figure  3-4). 

The  next  set  of  runs  where  designed  to  see  whether  the  estimated 
covariance  matrix  actually  reflects  the  observed  error.  The  expected  value  of 
the  observed  and  estimated  ECI  position  is  approximated  from  the  simulation 
to  be  the  observed  position  averaged  over  50  runs  and  the  estimated  position 
averaged  over  50  runs.  Additionally,  the  actual  covariance  matrix  of  the  error 
in  position  can  be  estimated  from  the  simulation  by  computing  the  empirical 
variance  and  covariance  of  the  errors  from  the  50  simulation  runs.  The 
average  of  the  50  estimated  covariance  matrices  will  be  an  unbiased  estimator. 

The  ultimate  goal  is  to  find  a  single  value  which  reflects  the  quality  of  a 
sensor  platform  observations.  Thus  the  total  variance  with  respect  to  ECI 
position  was  chosen  to  be  a  measure  of  quality  of  a  sensor  platform  observation, 
where  the  total  variance  is  defined  to  be: 

3  3  3 

V(U)  =  V  V(Yi)  +2x£  T  COV(YjYj) 

1*1  i-i  j-i 

where  Y,  i  =  1,2,3  is  the  ECI  X,  Y  and  Z  positional  variance,  respectively, 
and  U  =  Y,  +  Y2  +  Y3. 
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Figure  3-4  Normal  probability  plots  of  the  error  distributions. 
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The  ECI  position  covariance  matrix  was  chosen  as  opposed  to  the  launch 
parameter  covariance  because  the  launch  parameter  is  expressed  in  geographic 
coordinates,  time,  altitude,  azimuth  and  pitch  angle.  It  was  felt  that  it  would 
be  difficult  to  transform  this  data  into  a  common  variance.  Since  the  estimated 
total  variance  was  chosen  to  represent  the  algorithm’s  error,  there  should  a 
non-decreasing  relationship  between  the  observed  total  variance  in  position  and 
the  estimated  total  variance  in  position  derived  from  the  tracking  algorithm. 
The  error  in  position  (measured  in  radial  distance  from  the  estimated  and  real 
booster  position)  should  also  be  monotonically  related  to  observed  and 
estimated  total  variance. 

The  simulation  results,  observed  total  variance  (OTVAR)  vs.  expected  total 
variance  (ETVAR)  and  radial  error  (Figure  3-5  and  Figure  3-6)  were  plotted. 

The  estimated  total  variance  consistently  underestimates  the  observed 
total  variance.  Additionally,  between  approximately  time  20  and  time  50  the 
ETVAR  does  not  represent  OTVAR  or  the  radial  error.  The  average  radial 
error  is  fairly  well  estimated  by  OTVAR  (Figure  3-7). 

Additional  analysis  of  the  average  values  of  the  observed  and  estimated 
variances  indicated  that  the  values,  using  paired  comparisons,  were  not 
statistically  different.  The  estimated  covariance  was  of  the  order  of  10-6, 
significantly  smaller  in  value  than  the  observed  covariance  (order  of  101).  This 
indicates  that  the  tracking  algorithm  greatly  underestimates  the  covariance, 
and  is  an  explanation  as  to  why  ETVAR  underestimates  OTVAR. 


The  non-linear  relationship  of  ETVAR  and  OTVAR  and  the  radial  error 
over  time  indicates  that  there  was  a  failure  in  the  templates  to  accurately 
model  the  booster  trajectory  (template  mis-match)  or  that  the  data  set  for  the 
trajectory  was  faulty  (Figure  3-5). 

These  results  were  discussed  with  the  algorithms  author,  Nelson 
Rasmussen,  Martin  Marietta,  and  the  following  suggestions  were  made: 

1.  The  error  that  was  observed  during  the  simulation  run  time  from 
approximately  20  to  50  seconds  corresponds  to  missile  pitch  over  in  the 
booster’s  flight  profile.  It  might  be  the  case  that  the  templates  do  not 
model  this  well  causing  template  mis-match. 

2.  When  template  mismatching  occurs,  the  tracking  algorithm  might  not 
converge  well  and  will  produce  error. 

3.  An  experiment  was  run  in  which  the  initial  guess  of  the  launch  parameter 
was  varied  and  the  convergent  points  compared.  Out  of  12  different 
initial  guesses,  10  convergent  points  were  observed. 

The  experiment  described  by  Nelson  Rasmussen,  was  again  tested  on  the 
tracking  algorithm.  The  launch  azimuth  was  varied  from  0  to  2t  from  its 
initial  value  in  30  degree  increments  and  the  resulting  estimated  launch 
parameter  state  vectors  were  compared.  It  was  observed  during  a  large 
proportion  of  the  time  that  the  algorithm  converged  to  different  values,  but 
that  the  values  were  very  close  (the  same  out  to  8  decimal  places).  The 
geographic  positions  were  for  all  practical  purposes  the  same,  only  the 
launch  azimuth  and  pitch  parameters  converged  to  different  values. 

This  would  indicate  that  the  algorithm  is  robust  to  slight  changes  and 
errors  in  the  entering  arguments  or  initial  guess  of  the  launch  parameters. 
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ESTIMATED  VARIANCE  VS.  OBSERVED  VARIANCE 


30WIHVA  1V101  Q3Ad3S0O 
Figure  3-5  Plot  of  ETVAR  vs  OTVAR. 
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ETVAR,  OTVAR,  AND  AVERAGE  ERROR 


0 

H0HH3  30VH3AV  ONV  HVAIO 
Figure  3-8  Plot  of  ETVAR  vs  OTVAR  and  the  average  radial  error. 
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ETVAR 


OBSERVED  TVAR  VS  ERROR 


Since  the  algorithm  appears  to  work  for  most  time  periods,  it  was  believed 
that  some  template  mis-matching  was  occurring  during  the  time  period  from 
20  to  50  seconds. 

B.  SENSOR  IN  A  CONSTELLATION 

An  individual  sensor  platform  may  in  many  cases  not  give  good 
information.  Due  to  differences  in  viewing  angle,  sun  back-lighting  and  distance 
from  the  booster,  messages  from  different  sensors  v  ill  generate  messages  of 
varying  quality.  This  requires  analysis  of  the  message  quality  from  a 
constellation  of  sensors. 

If  a  population  of  sensor  platforms  observe  a  booster  launch  and  initiates 
a  track,  how  will  the  error  of  one  sensor  platform  compare  to  the  error  of 
another  sensor  tracking  the  same  booster? 

If  a  population  of  sensors  is  tracking  a  booster,  it  would  be  expected  that 
sensors  with  good  track  information  would  have  smaller  observed  radial  error, 
as  measured  by  the  estimated  ECI  position  covariance  matrix  (i.e.  smaller 
distance  between  the  estimated  position  and  the  real  position).  The  total  ECI 
(ETVAR)  is  again  used  as  a  single  measure  of  information. 

The  hypothesis  of  no  correlation  between  populations  of  ranked  pairs  of 
the  observed  error  and  the  estimated  total  variance  can  be  tested  by  using 
Spearman's  Rank  Correlation  Test  (Mendenhall,  1989).  The  rank  correlation 
coefficient  r,  is  calculated  by  using  the  ranks  of  the  paired  measurements  of  the 
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two  variables,  observed  error  and  estimated  variance.  Let  the  observed  error 
of  the  Ith  run  be  Xj  and  estimated  variance  of  the  i*  run  be  Y„  for  i  =  1,..,25  the 
correlation  r,  is  calculated  by: 

25 

eV  (Xi-Yj)2 

r  =  1--M - 

25  (252-l) 

Where  25  represents  the  number  of  sensors  tracking  the  booster.  This  assumes 


c  = 


v  =  n-2 


that  there  are  no  ties  in  either  the  X*  or  Ys  observations.  The  statistic  r,  can 
be  used  to  generate  a  t-statistic  used  for  hypothesis  testing  (Kendall,  1990). 

The  t-statistic  can  then  be  used  to  test  the  following  hypothesis  at  each 
time  period: 

1.  Null  Hypothesis:  Ha:  There  is  no  association  between  the  ranked  pairs. 

2.  Alternative  Hypothesis:  Ha:  The  correlation  between  the  ranked  pairs 
is  positive. 

A  simple  constellation  was  constructed  to  test  this  hypothesis.  To  reduce 
variation  introduced  into  the  simulation  testing,  the  constellation  was  centered 
on  the  ECI  position  of  the  sensor  used  in  all  pervious  experiments.  Each 
sensor  was  2  degrees  off  azimuth  of  neighbor  (There  are  cases  when  the 
algorithm  will  fail  to  generate  a  track.  To  minimize  these  errors,  it  was 
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decided  to  use  a  constellation  in  which  the  sensor  platforms  would  be  assured 
of  generating  a  track).  This  constellation  is  for  illustrative  purpose 
(Figure  3-8). 

The  tracking  algorithm  was  run  using  these  initial  positions.  The  estimated 
positions,  real  positions  and  covariance  matrix  output  data  was  run  through  a 
post-processor  (Appendix  2).  The  post  processor  calculated  the  radial  error 
(norm  of  estimated  position  and  real  position)  and  total  estimated  variance 
(ETVAR)  for  each  time  period.  These  variables  were  then  used  for  the  rank 
correlation  test.  The  results  indicate  no  correlation  between  the  observed  error 
and  the  estimated  total  variance.  In  most  cases  the  t-statistic  did  not  reject  the 
Null  Hypothesis  (Table  I).  At  the  95%  level  of  significance  the  t-critical  value 
with  23  degrees  of  freedom  is  2.069 

An  attempt  was  made  to  improve  the  correlation  between  the  radial  error 
and  the  variance  by  modeling  the  error  as  a  function  of  the  variance  covariance 
matrix.  If  the  error  can  be  predicted  by  the  covariance  matrix,  this  should 
improve  the  correlation  and  perhaps  give  an  absolute  measure  of  error. 

The  following  model  was  used: 

where  Y  is  the  radial  error  and  Xj  j  =  1,..,6  are  the  entries  from  the  covariance 
matrix. 
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Figure  3*7  Simple  constellation  of  25  SBls. 
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Least  Squares  regression  was  used  to  estimate  the  values  of  0.  The 
regression  will  only  show  what  are  important  factors  in  forecasting  error  and 
estimate  the  coefficients. 

When  using  the  ECI  covarinace  matrix,  the  coefficients  will  only  be 
optimal  for  a  limited  geographic  location  (unless  the  coefficient  are  the  same 
for  the  X,  Y,  and  Z,  variance  and  covariance).  For  this  reason  this  analysis  was 
additionally  tested  on  the  launch  parameter  variance.  These  coefficients  should 
be  able  to  be  used  globally. 

Three  cases  were  looked  at  for  the  regression  model.  The  launch  parameter 
variance  and  ECI  covariance  from  time  0  to  100,  and  the  ECI  covariance  from 
time  0  to  36  and  from  time  38  until  100  (Table  II).  The  time  periods  0  to  36 
and  38  to  100  were  chosen  to  observe  the  behavior  of  the  model  prior  to 
booster  pitch  over  and  during  the  ballistic  phase  of  the  flight. 

In  the  launch  parameter  covariance  model,  only  the  azimuth  and  pitch 
variance  were  used.  The  variance  for  latitude  and  longitude  were 
approximately  zero,  the  altitude  variance  was  constant,  and  the  use  of  time 
variance  resulted  in  even  a  lower  value  of  R. 

These  coefficients  were  then  used  in  the  post-processor  to  calculate  the 
value  of  ETVAR.  The  post  processor  then  computed  the  Spearman’s  Rank 
Correlation  Test,  as  previously  described  (see  Appendix  2). 

This  method  of  choosing  the  coefficient  to  calculate  ETVAR  did  not  change 
the  results:  in  most  cases  it  was  not  possible  to  reject  the  Null  Hypothesis 
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Table  I  RESULTS  OF  SPEARMAN'S  RANK  CORRELATION  TEST, 
ESTIMATED  VARIANCE  CALCULATED  FROM  ECI  COVARIANCE  MATRIX, 
&  “  1,1, 1,2, 2, 2 


TUB 

r. 

t - stat 

d.f . 

6 

0.05154 

0.2475 

23 

8 

0.32462 

1.6459 

23 

10 

0.03538 

0.1698 

23 

12 

•0.06462 

-0.3105 

23 

14 

0.00000 

0.0000 

23 

16 

-0.14693 

-0.7123 

23 

18 

0.04768 

0.2287 

23 

20 

0.16846 

0.8196 

23 

22 

0.38923 

2.0264 

23 

24 

-0.03154 

-0.1513 

23 

26 

0.03692 

0.1772 

23 

28 

0.30385 

1.5295 

23 

30 

-0.29231 

-1.4658 

23 

32 

-0.09692 

-0.4670 

23 

34 

0.05000 

0.2401 

23 

36 

0.18308 

0.8931 

23 

38 

-0.23846 

-1.1776 

23 

40 

-0.05154 

-0.2375 

23 

42 

0.06769 

0.3254 

23 

44 

-0.07692 

-0.3701 

23 

46 

0.00769 

0.0368 

23 

48 

•0.06154 

-0.2957 

23 

50 

-0.15154 

-0.7352 

23 

52 

0.25385 

1.2586 

23 

54 

0.03385 

0.1624 

23 

56 

-0.35846 

-1.8415 

23 

58 

-0.06000 

-2.8827 

23 

60 

-0.17769 

-0.8659 

23 

62 

0.21000 

1.0300 

23 

64 

-0.05769 

-0.2771 

23 

66 

-0.00462 

-0.0221 

23 

68 

0.08000 

0.3849 

23 

70 

-0.27000 

-1.3448 

23 

72 

0.24077 

1.1896 

23 

74 

0.39846 

2.0835 

23 

76 

0.01385 

0.0664 

23 

78 

0.26308 

1.3077 

23 

80 

-0.01538 

-0.0737 

23 

82 

-0.29308 

-1.4701 

23 

84 

-0.10468 

-0.5045 

23 

86 

0.25385 

1.2586 

23 

88 

-0.17077 

-0.8312 

23 

90 

-0.17000 

-0.8273 

23 

92 

0.10077 

0.4857 

23 

94 

-0.07846 

-0.3774 

23 

96 

-0.19846 

-0.9711 

23 

98 

0.06615 

0.3179 

23 

100 

-0.36769 

-1.8962 

23 

25 


1 


(Tables  II-VII).  Using  the  launch  parameter  covariance  to  forecast  error 
produced  higher  correlation  and  in  some  cases  the  Null  Hypothesis  could  be 
rejected. 

Thus,  there  are  cases  when  the  launch  parameter  covariance  matrix  could 
be  used  an  indicator  of  the  quality.  However,  positive  correlation  does  not 
guarantee  high  quality  messages  or  the  ability  to  forecast  error.  A  good  model 
that  forecast  error  would  result  in  a  rank  correlation  that  is  close  to  1.0,  and 
a  t-statistic  that  is  significantly  large. 
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Table  II  RESULTS  OF  SPEARMAN'S  RANK  CORRELATION  TEST, 
ESTIMATED  VARIANCE  CALCULATED  FROM  ECI  COVARIANCE  MATRIX, 
0  *  -390,172,-69,-1024,125200,-387470 


TKE 

r. 

t-stat 

d.f . 

6 

0.04231 

0.2030 

23 

8 

•0.30615 

•1.5423 

23 

10 

0.18154 

0.8853 

23 

12 

-0.16615 

-0.8088 

23 

14 

0.11615 

0.5608 

23 

16 

0.10692 

0.5157 

23 

18 

-0.16000 

-0.7774 

23 

20 

0.32000 

1.6198 

23 

22 

0.09077 

0.4371 

23 

24 

-0.18308 

-0.8931 

23 

26 

0.11692 

0.5646 

23 

28 

-0.01692 

-0.0811 

23 

30 

0.25846 

1.2831 

23 

32 

0.21462 

1.0538 

23 

34 

0.45769 

2.4687 

73 

36 

-0.06846 

-0.3291 

23 

38 

-0.14154 

-0.6857 

23 

40 

0.17924 

0.8737 

23 

42 

-0.01462 

-0.0701 

23 

44 

-0.29846 

-1.4997 

23 

46 

-0.23000 

-1.3343 

23 

48 

-0.19231 

-0.9398 

23 

50 

-0.00077 

-0.0036 

23 

52 

0.13923 

0.6743 

23 

54 

-0.20154 

-0.9867 

23 

56 

-0.31538 

-1.5939 

23 

58 

-0.09077 

-0.4371 

23 

60 

-0.25462 

-1.2627 

23 

62 

-0.07692 

-0.3700 

23 

64 

0.17846 

0.8698 

23 

66 

0.04462 

0.2141 

23 

68 

-0.22769 

-1.1214 

23 

70 

0.35462 

1.8188 

23 

72 

0.13385 

0.6477 

23 

74 

0.02308 

0.1107 

23 

76 

0.13769 

0.6667 

23 

78 

-0.06077 

-0.2919 

23 

80 

0.14077 

0.6819 

23 

82 

0.06077 

0.2919 

23 

84 

0.09231 

0.4445 

23 

86 

0.31538 

1.5938 

23 

88 

-0.02538 

-0.1217 

23 

90 

0.12462 

0.6023 

23 

92 

-0.14308 

-0.6933 

23 

94 

-0.18923 

-0.9242 

23 

96 

-0.23923 

-1.0181 

23 

98 

0.03000 

0.1439 

23 

100 

-0.36154 

-1.8597 

23 
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Tabl*  III  RESULTS  OP  SPEARMAN'S  RANK  CORRELATION  TEST, 
ESTIMATED  VARIANCE  CALCULATED  FROM  SCI  COVARIANCE  MATRIX, 
TIME  0  -  36,  fi  m  -308,9,-370,84,-1026,243 


TBB 

r. 

t-atat 

d.f . 

6 

0.29846 

1.49973 

22 

8 

0.30769 

1.55088 

22 

10 

0.20692 

1.01432 

23 

12 

0.26077 

1.29543 

23 

14 

0.38308 

1.98889 

23 

16 

0.25615 

1.27087 

23 

18 

0.27538 

1.37382 

23 

20 

0.26692 

1.32831 

23 

22 

0.34462 

1.76056 

23 

24 

0.40385 

2.11710 

23 

26 

0.25769 

1.27905 

23 

28 

0.48154 

2.63500 

23 

30 

0.24385 

1.20584 

23 

32 

0.16923 

0.82348 

23 

34 

0.34077 

1.73832 

23 

36 

0.13462 

0.65152 

23 

38 

0.07231 

0.34769 

23 

40 

0.01385 

0.06641 

23 

42 

0.14308 

0.69331 

23 

44 

-0.06385 

-0.30682 

23 

46 

0.37462 

1.93769 

23 

48 

0.34154 

1.74276 

23 

50 

0.41538 

2.18999 

23 

52 

0.29077 

1.45745 

23 

54 

0.32231 

1.63287 

23 

56 

0.11385 

0.54956 

23 

58 

0.26308 

1.30774 

23 

60 

0.11692 

0.56462 

23 

62 

0.15692 

0.76202 

23 

64 

0.28462 

1.42386 

23 

66 

0.19615 

0.95936 

23 

68 

0.23538 

1.16150 

23 

70 

0.29923 

1.50397 

23 

72 

0.08615 

0.41472 

23 

74 

0.05692 

0.27344 

23 

76 

0.10154 

0.48949 

23 

78 

0.44615 

2.39082 

23 

80 

-0.04769 

-0.22898 

23 

82 

0.55077 

3.16465 

23 

84 

-0.09077 

-0.43712 

23 

86 

0.00154 

0.00738 

23 

88 

0.20462 

1.00251 

23 

90 

0.10308 

0.49699 

23 

92 

0.18308 

0.89310 

23 

94 

0.10231 

0.49324 

23 

96 

0.42077 

2.22444 

23 

98 

-0.10769 

-0.51950 

23 

100 

0.05923 

0.28456 

23 

Table  IV  RESULTS  OF  SPEARMAN'S  RANK  CORRELATION  TEST, 
ESTIMATED  VARIANCE  CALCULATED  FROM  ECI  COVARIANCE  MATRIX 
TIME  38  -  100,  0  -  663,342,2537,-1031,-4036,2080 


TZB 

r. 

t  -  stat 

d.f . 

6 

0.12769 

0.61745 

22 

8 

0.39154 

2.04068 

22 

10 

0.17538 

0.85436 

23 

12 

0.14231 

0.68950 

23 

14 

0.15615 

0.75819 

23 

16 

0.26308 

1.30774 

23 

18 

0.28308 

1.41549 

23 

20 

0.33154 

1.68532 

23 

22 

0.44846 

2.40629 

23 

24 

-0.14769 

-0.71616 

23 

26 

0.23231 

1.14545 

23 

28 

0.16385 

0.79654 

23 

30 

-0.08308 

-0.39981 

23 

32 

0.09692 

0.46703 

23 

34 

0.30077 

1.51247 

23 

36 

0.14308 

0.69331 

23 

38 

0.48462 

2.65698 

23 

40 

0.19077 

0.93201 

23 

42 

0.59000 

3.50450 

23 

44 

0.43385 

2.30930 

23 

46 

0.70231 

4.73139 

23 

48 

0.35692 

1.83244 

23 

50 

0.53385 

3.02778 

23 

52 

0.26538 

1.32007 

23 

54 

-0.03538 

-0.16981 

23 

56 

0.61692 

3.75930 

23 

58 

0.47538 

2.59141 

23 

60 

0.28462 

1.42386 

23 

62 

-0.14462 

-0.70092 

23 

64 

-0.13538 

-0.65531 

23 

66 

-0.00692 

-0.03320 

23 

68 

0.16231 

0.78886 

23 

70 

0.24923 

1.23422 

23 

72 

-0.09615 

-0.46328 

23 

74 

0.45462 

2 .44784 

23 

76 

0.01615 

0.07748 

23 

78 

0.04538 

0.21788 

23 

80 

0.18846 

0.92032 

23 

82 

0.00769 

0.03689 

23 

84 

-0.03538 

-0.16981 

23 

86 

0.41000 

2.15582 

23 

88 

-0.05385 

-0.25861 

23 

90 

0.12692 

0.61366 

23 

92 

0.16923 

0.82346 

23 

94 

0.08846 

0.42592 

23 

96 

0.36308 

1.86878 

23 

98 

0.18462 

0.90087 

23 

100 

0.26385 

1.31185 

23 
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Tabl«  V  CORRELATION  OF  RADIAL  ERROR  TO  LAUNCH  PARAMETER 
VARIANCE ,  REGRESSION  ESTIMATES  OF  fi  -  0, 0, 0, 0f 1.3, -2 . 6 


TUB 

r. 

t-stat 

d.f 

6 

0.38615 

2.00765 

23 

8 

0.38769 

2.01706 

22 

10 

0.30308 

1.52524 

23 

12 

0.11385 

0.54956 

23 

14 

0.34846 

1.78291 

23 

16 

0.18154 

0.88534 

23 

18 

0.26077 

1.29543 

23 

20 

0.10846 

0.52325 

23 

22 

0.28077 

1.40296 

23 

24 

0.22538 

1.10945 

23 

26 

0.13769 

0.66670 

23 

28 

0.15769 

0.76585 

23 

30 

0.36231 

1.86422 

23 

32 

0.27615 

1.37797 

23 

34 

0.08385 

0.40353 

23 

36 

0.34000 

1.73388 

23 

38 

0.22385 

1.10148 

23 

40 

0.29692 

1.49125 

23 

42 

0.16846 

0.81963 

23 

44 

0.17615 

0.85822 

23 

46 

0.09462 

0.45580 

23 

48 

0.24231 

1.19776 

23 

50 

0.01923 

0.09224 

23 

52 

0.36846 

1.90082 

23 

54 

0.34154 

1.74276 

23 

56 

0.20308 

0.99465 

23 

58 

0.22692 

1.11744 

23 

60 

0.16846 

0.81963 

23 

62 

0.54154 

3.08933 

23 

64 

0.09615 

0.46328 

23 

66 

0.03154 

0.15133 

23 

68 

-0.06077 

•0.29198 

23 

70 

0.30154 

1.51673 

23 

72 

0.37154 

1.91922 

23 

74 

0.29462 

1.47855 

23 

76 

0.40077 

2.09787 

23 

78 

0.20308 

0.99465 

23 

80 

-0.10923 

-0.52701 

23 

82 

0.18462 

0.90087 

23 

84 

0.24308 

1.20180 

23 

86 

0.33923 

1.72944 

23 

88 

0.24000 

1.18565 

23 

90 

0.30154 

1.51673 

23 

92 

-0.05692 

-0.27344 

23 

94 

0.03154 

0.15133 

23 

96 

0.04462 

0.21418 

23 

98 

-0.06769 

-0.32539 

23 

100 

0.12231 

0.59100 

23 
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Table  VI  CORRELATION  OF  RADIAL  ERROR  TO  LAUNCH  PARAMETER 
VARIANCE,  ESTIMATED  COEFFICIENTS  FOR  TIME  0-36, 

0  -  0,0,0,0,-92.3,0.33 


7SB 

r. 

t-stat 

d.f 

6 

0.38615 

2.00765 

23 

8 

0.38769 

2.01706 

22 

10 

0.30308 

1.52524 

23 

12 

0.11385 

0.54956 

23 

14 

0.19308 

0.94372 

23 

16 

0.26462 

1.31596 

23 

18 

0.32385 

1.64158 

23 

20 

0.51846 

2.90779 

23 

22 

0.43923 

2.34477 

23 

24 

0.52385 

2.94933 

23 

26 

0.25385 

1.25863 

23 

28 

0.17846 

0.86984 

23 

30 

0.08538 

0.41099 

23 

32 

0.18231 

0.88922 

23 

34 

0.18769 

0.91643 

23 

36 

0.31154 

1.57234 

23 

38 

0.19538 

0.95545 

23 

40 

0.50231 

2.78595 

23 

42 

0.36231 

1.86422 

23 

44 

0.15231 

0.73906 

23 

46 

0.19923 

0.97502 

23 

48 

0.26462 

1.31596 

23 

50 

0.27462 

1.36967 

23 

52 

0.09846 

0.47451 

23 

54 

0.20923 

1.02615 

23 

56 

0.34615 

1.76949 

23 

58 

0.39769 

2.07872 

23 

60 

0.45846 

2.47403 

23 

62 

0.09231 

0.44459 

23 

64 

0.15615 

0.75819 

23 

66 

0.24308 

1.20180 

23 

68 

0.10923 

0.52701 

23 

70 

0.20923 

1.02615 

23 

72 

0.25923 

1.28723 

23 

74 

0.43769 

2.33460 

23 

76 

0.41000 

2.15582 

23 

78 

0.42769 

2.26915 

23 

80 

0.54538 

3.12052 

23 

82 

0.19538 

0.95545 

23 

84 

0.62923 

3.88267 

23 

86 

0.39923 

2.08828 

23 

88 

0.28923 

1.44903 

23 

90 

0.25154 

1.24641 

23 

92 

0.24154 

1.19372 

23 

94 

0.13615 

0.65911 

23 

96 

-0.05538 

-0.26602 

23 

98 

0.34000 

1.73388 

23 

100 

0.46077 

2.48983 

23 

31 


C.  USE  OF  MODEL 


In  practical  terms,  one  would  like  to  prune  messages  of  poor  quality  from 
the  constellation.  Pruning  messages  could  greatly  reduce  the  time  required  for 
the  remaining  messages  to  propagate  through  the  constellation.  The  process 
of  message  propagation  will  be  defined  as  flooding. 

An  experiment  was  conducted  to  estimate  the  distribution  of  the  quality 
of  messages  generated  in  the  test  constellation.  The  model  that  consistently 
produced  the  highest  correlation  and  t-statistic  was  used  to  estimate  the  quality 
of  the  messages:  Q  =  -92.33(Var(azimuth))  +  0.33(Var(pitch)). 

A  level  of  80  percent  pruning  was  simulated.  For  each  message  period, 
only  the  messages  with  the  5  lowest  values  of  Q  were  retained.  When  one  of 
the  top  5  messages  was  present  among  the  pruned  messages  a  value  of  1  was 
given,  otherwise,  0.  Twenty  independent  runs  were  conducted. 

The  resulting  matrixes  were  summed  together  and  the  entries  were 
divided  by  20.  This  gives  the  proportion  of  times  that  one  of  the  five  highest 
quality  messages  was  present  at  a  given  quality  level  (Figure  3-9,  3-10).  If  Q 
forecast  error  well,  the  proportions  for  zero  through  five  of  Figures  3-9  and  3-10 
would  be  1.0,  and  the  proportions  from  six  to  twenty  five  would  be  0.0. 

The  proportions  for  the  five  highest  quality  messages  present  in  the 
primed  messages  are  less  than  one,  indicating  that  lower  quality  messages 
would  be  present  in  the  pruned  messages.  This  shows  that  Q  does  not 
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Figure  3-9  Density  plot  of  the  proportion  of  times  that  the  5  best 


messages  are  present  at  an  estimated  quality. 
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estimate  the  error  well  and  suggest  that  use  of  the  estimated  launch  parameter 
covariance  matrix  is  not  a  good  decision  rule  for  message  pruning. 

The  proportion  of  high  quality  messages  that  are  present  is  not  high  or 
constant,  as  indicated  in  Figures  3-9  and  Figure  3-10.  This  indicates  that  using 
the  launch  covariance  matrix  to  estimate  quality  does  not  significantly  improve 
the  chance  of  identifying  a  high  quality  message.  In  many  cases,  it  is 
detrimental.  From  scan  20  to  scan  35,  corresponding  to  time  46  to  76,  the 
estimated  quality  is  worse  than  if  messages  are  pruned  randomly.  This  is 
similarly  reflected  in  the  results  of  the  rank  correlation  test  (Table  XI). 

In  the  cases  tested,  the  launch  parameter  covariance  or  the  ECI 
covariance  did  not  well  represent  the  observed  radial  error.  This  makes  it 
difficult  to  constantly  forecast  radial  error  or  use  this  information  to  make  good 
decisions  in  message  pruning. 
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Figure  3*10  Contour  plot  of  the  proportion  of  5  best  messages  are  present 
at  a  quality  level 
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IV.  AZIMUTH  AND  ELEVATION  VARIANCE 

A.  AZIMUTH  AND  ELEVATION 

The  error  generated  in  the  tracking  algorithm  can  be  localized  to  two 
general  causes:  1.  template  mis- matching,  2.  variance  in  the  azimuth  and 
elevation.  By  reducing  the  variance  in  the  azimuth  and  elevation  measure  one 
would  expect  the  radial  error  and  the  variance  in  the  radial  error  to  be  reduced. 

Additionally,  by  reducing  the  azimuth  and  elevation  error,  the  amount  of 
noise  in  the  system  might  be  reduced  to  the  point  where  it  would  be  able  to 
forecast  the  error  from  the  estimated  ECI  covariance  matrix. 

The  variance  of  the  sensor  azimuth  and  elevation  was  changed  from  a 
value  we  will  call  A  to  .33A  and  .16A.  The  tracking  algorithm  was  then  run 
using  the  constellation  of  25  sensors  observing  a  single  booster  over  the  time 
from  0  to  100  seconds. 

Spearman’s  Rank  Correlation  Test  was  used  to  test  the  hypothesis  of 
correlation  between  observed  radial  error  and  the  estimated  total  variance. 
The  analysis  was  similar  in  design  to  that  conducted  on  previous  runs  using  the 
azimuth  and  elevation  variance  of  A.  The  regression  model  was  used  to  find 
the  estimates  of  the  coefficients  that  best  forecast  the  error. 

For  both  the  experiment  using  the  variance  set  at  .33A  and  .16A,  it  was 
not  possible  to  reject  the  Null  Hypothesis  (Table  VIII,  IX).  It  must  be  noted 
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that  at  certain  times  (after  time  80)  the  algorithm  failed  to  report  a  position. 
This  may  be  caused  by  the  information  matrix  becoming  singular,  making  it 
impossible  to  invert.  This  problem  was  more  pronounced  when  the  smaller 
variance  value  were  used. 

B.  MEAN  RADIAL  ERROR 

Reducing  the  variance  in  the  azimuth  and  elevation,  all  other  factors  being 
the  same,  should  reduce  the  tracking  system  error.  This  would  result  in 
smaller  radial  error  of  the  tracking  system  and  may  have  implication  in 
choosing  specification  for  sensor  performance. 

The  mean  radial  error  and  standard  deviation  were  calculated  for  the 
three  cases:  variance  in  the  azimuth  and  elevation  set  at  A,  .33A,  and  .16A. 
All  three  cases  used  25  sensors  observing  a  single  booster.  The  95  percent 
confidence  interval  for  each  time  period  was  then  calculated  and  plotted  (Figure 
4-9). 

The  graphical  results  show  that  in  most  cases  there  is  a  significant 
difference  in  the  mean  radial  error  at  different  levels  of  azimuth  and  elevation 
variance.  Additionally,  the 
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Table  IX  CORRELATION  OF  RADIAL  ERROR  TO  ESTIMATED 
VARIANCE,  AZIMUTH  AND  ELEVATION  VARIANCE  =  .33A, 
0  =  180,-47,-6668,-4036,-136680,-347470 


TIME 

r. 

t-atat 

d.f. 

6 

0.18154 

0.88534 

22 

8 

0.18692 

0.91254 

23 

10 

-0.06692 

-0.32167 

23 

12 

0.31692 

1.60252 

23 

14 

-0.12385 

-0.59855 

23 

16 

-0.21692 

-1.06570 

23 

18 

•0.08308 

-0.39981 

23 

20 

0.20000 

0.97895 

23 

22 

-0.16692 

-0.81193 

23 

24 

0.24769 

1.22610 

23 

26 

0.16615 

0.80808 

23 

28 

-0.00462 

-0.02214 

23 

30 

-0.02769 

•0.13286 

23 

32 

0.37769 

1.95625 

23 

34 

0.26692 

1.32831 

23 

36 

0.04154 

0.19938 

23 

38 

-0.06462 

-0.31053 

23 

40 

0.05692 

0.27344 

23 

42 

-0.20385 

-0.99858  . 

23 

44 

0.01231 

0.05903 

23 

46 

0.10077 

0.48574 

23 

48 

0.26385 

1.31185 

23 

50 

-0.18769 

-0.91643 

23 

52 

-0.11231 

-0.54204 

23 

54 

0.02308 

0.11070 

23 

56 

-0.17231 

-0.83891 

23 

58 

-0.21308 

-1.04590 

23 

60 

0.01000 

0.04796 

23 

62 

0.01769 

0.08486 

23 

64 

0.29231 

1.46588 

22 

66 

0.20615 

1.01038 

22 

68 

0.27923 

1.39462 

22 

70 

0.05385 

0.25861 

22 

72 

0.21846 

1.07364 

22 

74 

0.19846 

0.97110 

22 

76 

0.28538 

1.42804 

22 

78 

0.24000 

1.18565 

22 

80 

0.73231 

5.15736 

18 

82 

0.86923 

8.43169 

16 

84 

model 

failure 

-  - 

86 

model 

failure 

-  - 

88 

model 

failure 

90 

model 

failure 

-  - 

92 

model 

failure 

-  - 

94 

model 

failure 

-  - 

96 

0.85769 

8.00008 

13 

98 

0.41077 

2.16068 

18 

100 

model 

failure 

-- 
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Table  X  CORRELATION  OF  RADIAL  ERROR  TO  ESTIMATE 
D  VARIANCE,  AZIMUTH  AND  ELEVATION  VARIANCE  =  .16 
A, 

0  =  -193,110,-3231,-2970,179640,-643690 


7SB 

r. 

t-stat 

6 

-0.04769 

•0.22898 

23 

8 

-0.08692 

-0.41845 

23 

10 

0.06154 

0.29569 

23 

12 

0.25385 

1.25863 

23 

14 

0.46615 

2.52694 

23 

16 

-0.06615 

-0.31796 

23 

18 

-0.09000 

-0.43338 

23 

20 

0.18000 

0.87758 

23 

22 

0.13154 

0.63637 

23 

24 

0.36077 

1.85512 

23 

26 

0.47846 

2.61314 

23 

28 

0.06000 

0.28827 

23 

30 

-0.18846 

-0.92032 

23 

32 

0.07462 

0.35884 

23 

34 

0.28692 

1.43643 

23 

36 

0.45692 

2.46353 

23 

38 

0.20385 

0.99858 

22 

40 

0.18615 

0.90864 

22 

42 

-0.27385 

-1.36552 

22 

44 

0.06769 

0.32539 

22 

46 

0.28385 

1.41967 

22 

48 

0.30923 

1.55945 

22 

50 

0.21615 

1.06174 

22 

52 

0.24538 

1.21394 

22 

54 

0.26154 

1.29953 

22 

56 

0.64846 

4.08529 

19 

58 

0.55846 

3.22868 

18 

60 

0.74231 

5.31294 

15 

62 

0.86538 

8.28212 

12 

64 

model 

failure 

-  - 

66 

model 

failure 

68 

model 

failure 

-  - 

70 

model 

failure 

72 

model 

failure 

-  - 

74 

model 

failure 

76 

model 

failure 

-  - 

78 

model 

failure 

80 

model 

failure 

82 

model 

failure 

-  - 

84 

model 

failure 

86 

model 

failure 

88 

model 

failure 

90 

model 

failure 

-  - 

92 

model 

failure 

94 

model 

failure 

96 

model 

failure 

-  . 

98 

model 

failure 

100 

model 

failure 

-  - 
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variance  in  the  observed  radial  error  increases  with  the  mean  radial  error. 
Note  that  if  data  was  not  available  (failure  to  report  a  position)  the  error  was 
reported  as  a  mean  of  zero  and  a  standard  deviation  of  zero. 

In  practical  terms,  the  more  accurately  a  measurement  of  azimuth  and 
elevation  can  be  made,  the  smaller  the  radial  error.  The  population  of 
messages  that  result  from  a  constellation  of  sensor  having  more  accurate 
measurements  will  have  messages  of  less  variance  and  of  more  consistent 
quality.  In  this  situation,  random  pruning  will  be  more  effective:  the  overall 
message  quality  is  higher,  and  at  any  given  pruning  level,  the  chances  are  that 
the  messages  will  have  similar  error. 
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95  PERCENT  CONFIDENCE  INTERVAL 


Figure  4-9  95  percent  confidence  interval  for  radial  error  at  Var  =  A, 
0.33A  and  0.16A. 
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TIME 


V.  CONCLUSIONS  AND  RECOMMENDATIONS 
The  Template  Based  Tracking  Algorithm  is  capable  of  estimating  the 
position  of  ballistic  bodies  or  boosters  with  just  a  single  optical  sensor.  An 
individual  sensor,  if  not  obscured  by  the  earth,  will  initiate  and  track  a 
thrusting  body  with  a  remarkable  degree  of  accuracy. 

However,  the  system  in  which  the  sensor  is  deployed  will  require 
numerous  sensors  in  a  constellation  orbiting  the  earth.  Any  launch  of  a  booster 
or  thrusting  body  will  be  viewed  by  a  number  of  sensors,  producing  a 
population  of  launch  parameter  messages  of  varying  degrees  of  quality. 

The  measure  of  quality  of  measurement  was  taken  to  be  the  radial  error, 
or  the  distance  from  the  estimated  position  to  the  real  position.  The  estimator 
of  message  quality  used  was  either  the  total  estimated  variance  calculated  from 
the  variance  of  the  sum  of  ECI  position  or  the  total  variance  calculated  when 
using  the  coefficients  estimated  by  using  least  square  regression  for  both  the 
ECI  covariance  matrix  and  the  launch  parameter  variance. 

The  correlation  between  the  observed  error  and  the  estimator  of  message 
quality  was  tested  using  Spearman’s  Rank  Correlation  Test.  It  was 
hypothesized  that  if  numerous  messages  are  generated  regarding  a  booster, 
messages  with  the  smallest  estimated  total  variance  will  have  the  best  quality 
launch  parameter  information  (smallest  radial  error). 
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Spearman’s  Rank  Correlation  Test  was  used  to  test  the  hypothesis  that 
radial  error  and  estimated  total  variance  were  positively  correlated.  In  most 
every  case  tested,  the  Null  hypothesis  could  not  be  rejected. 

Because  the  information  in  the  ECI  position  covariance  matrix  and  the 
launch  parameter  covariance  matrix  did  not  represent  the  observed  radial 
error,  a  deterministic  method  to  estimate  error  should  be  examined.  If  it  can 
be  determined  what  is  the  best  relative  viewing  angle  of  a  booster  from  a  given 
sensor  platform,  an  algorithm  could  be  generated  that  exploits  this  information. 

With  such  system,  the  launch  parameters  state  vector  from  those  sensor 
platforms  that  were  determined  to  have  the  smallest  error  radial  would  be 
allowed  to  flood  the  constellation.  Information  from  sensor  platforms  that  were 
in  a  poor  position  to  generate  booster  tracks  would  not  be  transmitted.  If 
practical,  this  method  could  reduce  the  queuing  problem  by  allowing  only 
information  of  high  quality  (determined  by  relative  viewing  angle  of  the  booster 
from  a  given  sensor)  to  be  transmitted  and  flooded  through  the  constellation. 

There  exists  a  problem  in  the  algorithm  during  the  time  period  from 
approximately  20  to  40  seconds.  This  failing  of  the  tracking  algorithm  to 
successfully  estimate  the  quality  of  a  message  could  be  the  result  of  several 
things.  After  discussion  with  Nelson  Rasmussen,  it  was  noted  that  the  a  priori 
information  derived  from  the  templates  might  not  model  that  portion  of  the 
flight  envelope.  This  would  cause  template  mis-match  and  result  in  a  degraded 
estimate  of  the  position.  A  final  point  might  be  that  the  time  increment  for 
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the  templates  is  too  coarse.  Thus  the  cubic  splines  would  not  represent  the 
changes  in  position,  velocity  and  acceleration  adequately. 

The  templates  have  time  increments  of  10  seconds  and  a  precision  of  2 
decimal  places.  The  actual  flight  data  of  the  boosters  has  3  decimal  places  of 
precision.  Because  the  acceleration,  altitude  and  downrange  are  so  variable 
during  the  initial  phase  of  flight,  it  may  be  possible  that  the  cubic  splines 
generated  from  the  templates  are  not  accurately  modeling  the  boosters 
trajectory.  It  is  proposed  that  the  templates  be  modeled  in  time  increments  of 
2  seconds  during  the  early  phase  of  flight  and  that  the  precision  be  at  least  3 
decimal  places.  This  alone  might  improve  the  template  fit,  resulting  in  better 
estimation  of  the  launch  parameters  and  a  corresponding  reduction  in  radial 
error. 

The  Template  Based  Tracking  Algorithm  operates  very  well.  Because  the 
algorithm  can  converge  quickly  to  a  sharp  answer,  it  is  felt  that  the  normally 
unimportant  second  order  effects  would  become  significant.  The  precision  and 
resolution  of  the  templates  may  contribute  significantly  to  the  observed  error. 
Additionally,  the  use  of  a  second  order  Taylor  expansion  to  estimate  Z  (azimuth 
and  elevation  measurement)  could  greatly  improve  the  ability  to  forecast  error. 

The  Template  Based  Tracking  Algorithm  procedure  to  solve  for  the  launch 
parameters  operates  well.  It  has  shown  the  ability  to  track  thrusting  bodies 
using  a  single  optical  sensor.  However,  at  this  time,  there  appears  to  be  little 
relationship  between  either  the  ECI  covariance  matrix  or  the  launch  parameter 


variance  and  the  radial  error.  Additional  testing  is  required  to  resolve  the 
problem  of  the  algorithm  ability  to  track  the  body  during  the  early  time  period 
associated  with  pitch-over,  from  approximately  20  to  60  seconds.  Once  this 
aspect  of  the  tracking  algorithm  is  adjusted,  it  may  be  possible  to  consistently 
forecast  radial  error.  At  the  present,  determining  quality  of  error  may  be 
likened  to  a  coin  toss. 
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APPENDIX  1 


PROGRAM  TVAR 


*  This  program  reads  in  the  real  eci  position,  estimated 

*  position  and  covariance  matrix  and  calculates  the  expected 

*  radial  error,  expected  value  of  the  observed  total  variance 

*  and  the  expected  value  of  the  estimated  total  variance  at 

*  each  time  interval  (scan).  It  takes  as  input,  50  runs 

*  from  a  single  sensor  observing  a  booster. 

*  by  Eric  Bechhoefer,  SMC  1089  NPS ,  Monterey  Ca  93943 

* 

IMPLICIT  NONE 


REAL* 8  ESTP0S(3 , 2500) 
REAL*8  REALPOS(3 , 2500) 
REAL* 8  EDATA(3 , 2500) 
REAL*8  RDATA(3 , 2500) 
REAL* 8  C0V(3 , 3 , 2500) 
REAL* 8  COW (3, 100) 

REAL* 8  C0VMTX(3 , 3 , 2500) 
REAL* 8  ECOVMTX(3 , 3 , 2500) 
REAL* 8  VAR (3 , 100) 

REAL* 8  EVAR(3 , 100) 

REAL* 8  EY1Y2(3 , 100) 

REAL* 8  MEAN (3, 100) 

REAL*8  SUMSQ(3 , 100) 

REAL* 8  SUM(3 , 100) 

REAL* 8  TOTV(IOO) 

REAL* 8  ETOTV(IOO) 
real*8  DIST(100) 

INTEGER  I.J.K.L.M.N 
INTEGER  START 
INTEGER  MARKER 
INTEGER  COUNT 
INTEGER  TIME 
INTEGER  DTIME 
PARAMETER  (DTIME  -  2) 
INTEGER  SLICES 


!  estimate  position  data 
!real  position  data 
! estimate  position  data 
!real  position  data 
!  covariance  data 


{calculate  covariance  data 
{observed  variance  data 
{estimated  variance  data 
{for  calculating  covariance 
{mean  distance  error 
{sum  square  error  data 
{sum  data  of  positional 
!  error 

! total  observed  variance 
{estimated  total  variance 
{radial  error 
! counters 

{start  of  an  array 
{marker 
!a  counter 
!sim  time 
{scan  rate 

•t 

{number  of  observations 


******begin 


OPEN  (20, FILE  -  'EST_ECI_P0S  \  STATUS  -  'OLD') 
OPEN  (25 .FILE  -  ' REAL_ECI_P0S  ',  STATUS  -  'OLD') 
OPEN  (50, FILE  -  'ECIJSOV  ',  STATUS  -  'OLD') 


OPEN  (UNIT  -  99,  FILE  -  'VARS  STATUS  -  'NEW') 

N  -  0 

print*, ' input  count,  which  is  the  number  or  repetitions' 
read* , count 

print*, ' input  slices  -' 
read*,  slices 
print*, ' input  start  time  ' 
read*, time 


TIME  -  6  +  time 
MARKER  -  COUNT  *  SLICES 
DO  10  I  -  1,  SLICES 
TOTV(I)  -  0.0 
ETOTV(I)  -  0.0 
dist(i)-0.0 
DO  5  J-l , 3 

SUM(J , I)  -  0.0 

SUMSQ(J.I)  -  0.0 
MEAN(J.I)  -  0.0 
EY1Y2(J , I)  -  0.0 
COW(J,I)  -  0.0 
EVAR(J , I)  -  0.0 

DO  3  M-1,3 

ECOVMTX(M, J , I+N)  -  0.0 
ECOVMTX (M , J , I+N+l )  -  0.0 
EC0VMTX(M, J , I+N+2)  -  0.0 
3  CONTINUE 

N  -  N  +  2 
5  CONTINUE 

10  CONTINUE 

DO  14  I  -  1 .MARKER 
DO  13  J  -  1,3 

RDATA(J.I)  -  0.0 
EDATA(J.I)  -  0.0 
DO  12  K  -  1,3 

COV(J ,K, I)  -  0.0 

12  CONTINUE 

13  CONTINUE 

14  CONTINUE 


CALL  FIX (COUNT , SLICES , RDATA , EDATA , COV) 


**finds  missing  data  points  in  the  input 


PRINT*, 'COMPLETED  READING  IN  ESTPOS,  REALPOS ,  AND  COV' 
DO  40  I  -  1, SLICES 
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DO  35  J  -  0,  COUNT-1 

L-( COUNT  *  I)-  (COUNT  -  1)  +  J 

DO  30  K  -  1,3 

ESTPOS(K, L)— EDATA(K, (I+(J*SLICES) ) ) 
REALPOS(K,L)-RDATA(K, (I+(J*SLICES) ) ) 

ESTPOS(K.L)  -  ESTPOS (K,L)  -  REALPOS(K.L) 

DO  25  M  -  1,3 

COVMTX(K,M, L)  -  C0V(K,M, (I+(J*SLICES) ) ) 

25  CONTINUE 

30  CONTINUE 

35  CONTINUE 
40  CONTINUE 

DO  70  I  -  0,  SLICES  -  1 

♦♦calculate  radial  error  and  mean  distance  error************ 

**and  expected  value  of  the  estimated  variance************** 

DO  60  J  -  1,  COUNT 

L  -  J  +  COUNT  *  I 
DO  50  K-1,3 

DIST(1+1)  -  ESTP0S(k, 1)**2  +  DIST(I+1) 

SUM(K, 1+1)  -  ESTPOS (K,L)  +  SUM(K,I+1) 

DO  45  M  -  1,3 

IF(K  .EQ.  M)THEN 

ECOVMTX(K,M, 1+1)  -  COVMTX(K,M,L)**2  + 

+  ECOVMTX(K,M, 1+1) 

ELSE 

ECOVMTX(K,M,I*l)  -  (COVMTX(K,M,L))  + 

+  ECOVMTX(K,M, 1+1) 

END  IF 

45  CONTINUE 

50  CONTINUE 

60  CONTINUE 
70  CONTINUE 

**f ind  mean  radial  error************************************ 

DO  90  I  -  1,  SLICES 

DIST(I)  -  SQRT(DIST(I)  /  REAL (COUNT) ) 

DO  80  J  -  1,  3 

MEAN(J , I)  -  SUM(J , I)  /  REAL(COUNT) 

DO  75  K  -  1,3 

ECOVMTX(J ,K, I)  -  EC0VMTX( J ,K, I)  / 

+  REAL (COUNT) 

75  CONTINUE 

80  CONTINUE 
90  CONTINUE 

DO  120  I  -  0,  SLICES  -  1 
DO  110  J  -  1,  COUNT 
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LrJ+50*I 
DO  100  K-1,3 

SUHSQ(K, 1+1)  -  ( ESTPOS (K , L)— MEAN (K , 1+1 ) ) **2 
+  +  SUMSQ(K, 1+1) 

100  CONTINUE 

110  CONTINUE 
120  CONTINUE 


DO  140  I  -  1,  SLICES 
DO  130  J  -  1,3 

VAR(J.I)  -  SUMSQ(J.I)  /  REAL (COUNT  -  1) 
130  CONTINUE 
140  CONTINUE 


**find  the  observed 


DO 


+ 

+ 

+ 

+ 

+ 

+ 


160  I  -  0,  SLICES  -  1 
DO  150  J  -  1,  COUNT 
L-J+50*I 

EY1Y2(1 , 1+1)  -  ( ESTPOS (l.L) -MEAN (1,1+1))  * 
(ESTPOS(2 ,L)-MEAN(2 , 1+1) )  + 
EY1Y2( 1,1+1) 

EY1Y2 ( 2 , 1+1 )  -  ( ESTPOS ( 2, L) -MEAN (2,1+1))  * 
( ESTPOS ( 2, L) -MEAN (2,1+1))  + 
EY1Y2(  2,1+1) 

EY1Y2(3, 1+1)  -  ( ESTPOS ( 1,L) -MEAN ( 3, 1+1))  * 
( ESTPOS ( 2 , L) — MEAN ( 2 , 1+1 ) )  + 
EY1Y2( 3,1+1) 


150  CONTINUE 
160  CONTINUE 

DO  180  I  -  1,  SLICES 
DO  170  J  -  1,  3 

COW(J,I)  -  EY1Y2(J ,  I)  /  REAL  (COUNT) 
170  CONTINUE 
180  CONTINUE 


**  FIND  TOTAL  VARIANCE  AS  THE  SUM  OF  ECI  X.Y.Z,  V(U)  -  *** 
**  V(X)+V(Y)+V(Z)  **************************************** 
DO  190  I  -  1,  SLICES 

TOTV(I)  -  VAR(l.I)  +  VAR(2 , I)  +  VAR(3,I)  + 

+  2*(COW(l ,  I)  +  COW  (2,1)  +  COW(3,I)) 

ETOTV(I)  -  ECOVMTX(l , 1,1)  +  ECOVMTX(2 , 2 , I)  + 

+  ECOVMTX(3 , 3,1)  +  2*(EC0VMTX(1 , 2 , I)  + 

+  ECOVMTX(l , 3,1)  +  ECOVMTX(2 ,3,1)) 

WRITE  (99,77)  TOTV(I) ,  ETOTV(I) ,dist(i) 
c  PRINT*,  TOTV(I),  ETOTV(I) 

77  FORMAT  (F13.4.4X,  F13 .4 ,4x, fl3 .4) 

190  CONTINUE 


STOP 

END 


APPENDIX  2 


PROGRAM  CROSS 

*  THIS  PROGRAM  WILL  READ  IN  THE  ESTIMATED  POSITION,  REAL 

*  POSITION  AND  THE  VARIANCE  COVARIANCE  MATRIX  FROM  THE 

*  SIMULATOR.  IT  WILL  THEN  FIND  MISSING  DATA,  TAKE  THE  ORDERED  *  DATA  AND 
CALCULATE  THE  RADIAL  ERROR  AND  THE  TOTAL  VARIANCE,  *  USE  SPEARMANS  RAND 
CORRELATION  test,  calculate  the 

*  t-  statistic  associated  with  that  correlation  .  by  Eric  R  *  Bechhoefer, 
SMC  1089,  NPS  Monterey,  Ca  93943 

* 


IMPLICIT  NONE 

************************************************************* 


REAL  EDATA(3 , 2500)  !  estimated  position  data 

REAL  RDATA(3 , 2500)  !  real  position  data 

REAL  COV(3 , 3 , 2500)  !  variance  data 

REAL  EVAR(100)  !  total  variance  from 


!  variance  matrix 


REAL  T_STAT 
REAL  CORR 

REAL  SUM,  SUMSQ,  MEAN,  STD 

REAL  ORDER_X(100) 

REAL  ORDER_Y ( 100 ) 

REAL  DIST(100) 

REAL  B1,B2,B3,B4,B5,B6 

INTEGER  I , J ,K, L,m,n,MM,MRK 

INTEGER  STIME 

INTEGER  DTIME 

INTEGER  START 

INTEGER  MARKER 

INTEGER  COUNT 

PARAMETER  (COUNT  -  25) 

INTEGER  SPACING 

INTEGER  SLICES 

INTEGER  BIN0(48 ,25) 

!  quality 


!  t-statistic 
!  correlation,  r, 

!  variable  for  calculating 
!  mean  radial  error  and  std. 
!  rank  order  of  radial  error 
!  rank  order  of  est.  var 
!  radial  error 
!  coefficients 
!  counters 
!  simulation  time 
!  scan  rate 
!  start  of  array 

!  number  of  sensor 

!  number  of  time  periods 
!  number  of  time  periods 
!  registers  a  1  if  a  top  5 
message  is  present 


character*8  a,b,c,d 
character*50  comment 

************************************************************ 
PRINT* ,' INPUT  THE  ESTPOS ,  REALPOS  AND  COV  MATRIX  NAME' 
READ* , A , B , C 

PRINT*, 'WHAT  IS  THE  NAME  OF  THE  OUTPUT' 

READ*,D 

PRINT*, 'ADD  ANY  COMMENTS?' 
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READ*,  COMMENT 


OPEN  (20, FILE  -  A  ,  STATUS  -  'OLD') 

OPEN  (25, FILE  -  B  ,  STATUS  -  'OLD') 

OPEN  (50, FILE  -  C  ,  STATUS  -  'OLD') 

OPEN  (60, FILE  -  D  ,  STATUS  -  'NEW') 

OPEN  (71, FILE  -  'DIST' .STATUS  -  'NEW') 

SLICES  -  48 
STIME  -  0 
SPACING  -  SLICES 

PRINT*, 'READ  IN  B1,B2,B3,B4,B5  AND  B6' 

READ* , B1 , B2 , B3 , £4 , B5 , B6 

PRINT*, 'IF  LAUGH  COVARIANCE  FILE,  M  -  1,  ELSE  ZERO' 
READ*, MM 

WRITE(60,*)  COMMENT 

WRITE (60,*) 'COEFFICIENTS  OF  B1,B2,B3,B4,B5,B6  ARE' 
WRITE(60,*)'B1-  ' ,B1, 'B2-' ,B2 
WRITE ( 60, *)'B3-  ',B3,'B4-  ',B4 
WRITE ( 60, *) 'B5-  ',B5,'B6-  ' ,B6 
62  FORMAT ( IX, 6F8.1) 


***  INITIALIZE  VARIABLES*******’ 
START  -  1 
STIME  -  STIME  +  6 
DTIME  -  2 
T_STAT  -  0.0 
MARKER  -  COUNT  *  SLICES 
DO  5  I  -  1, SLICES 
DO  3  J  -  1, COUNT 
BINO(I , J)  -  0 

3  CONTINUE 

5  CONTINUE 

DO  10  I  -  1,  COUNT 

ORDER_X(I)  -  REAL(I) 
ORDER  Y(I)  -  REAL(I) 
DIST(l)  -  0.0 
EVAR(I)  -  0.0 
10  CONTINUE 

DO  14  I  -  1, MARKER 
DO  13  J  -  1,3 

RDATA(J.I)  -  0.0 
EDATA(J.I)  -  0.0 
DO  12  K  -  1,3 

COV ( J , K , I )  -  0.0 

12  CONTINUE 

13  CONTINUE 

14  CONTINUE 

*****  FIXL  find  missing  data  **■* 
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CALL  FIX1( COUNT , SLICES , RDATA , EDATA , COV) 

PRINT*, 'COMPLETED  READING  IN  ESTPOS,  REALPOS,  AND  COV' 
*****  reorder  the  radial  error  and  total  variance  data*** 

DO  60  I  -  1, SPACING 
DO  50  J  -  0,  COUNT-1 
L  -  I  +  J*SLICES 
DO  40  K  -  1,3 

DIST(J+1)  -  DIST(J+1)+(EDATA(K,L) 

+  -RDATA(K,L))**2 

40  CONTINUE 

IF(MM  .NE.  1)THEN 

EVAR(J+1)  -Bl*COV ( 1 , 1 , L) **2  +B2*COV(2,2,L)**2 
+  B3*C0V ( 3 , 3 , L) **2  +B4(1)*C0V(1 , 2 ,L)  + 

B5 ( 1 ) *C0V ( 1 , 3 , L)  +  B6 ( 1 ) *C0V ( 2 , 3 , L) 

ELSE 

EVAR(J+1)  -  B1*C0V ( 1 , 1 , L) **2  +B2*COV(2 , 2 ,L)**2 
+B3*COV ( 3 , 3 , L) **2  +  B4(l)*COV(l,2,L)**2  + 
B5 ( 1 ) *COV ( 1 , 3 , L) **2  +  COV(2,3,L)**2 

ENDIF 

DIST(J+1)  -  SQRT(DIST(J+1) ) 

CONTINUE 

MRK  -  0 
SUM  -  0.0 
SUMSQ  -  0.0 
MEAN  -  0.0 

****  missing  data  is  identified,  reduces  d.f.  *************** 

DO  56  J-l, COUNT 

IF(DIST(J)  .It.  17 .0)THEN 

MRK  -  MRK  +  1 
DIST(MRK)  -  DIST(J) 

EVAR(MRK)  -  EVAR(J) 

SUM  -  SUM  +  DIST(J) 

ENDIF 

56  CONTINUE 

****  calculate  mean  radial  error  and  std.  ******************* 

IF  (MRK  .GT.  0)  MEAN  -  SUM/REAL(MRK) 

DO  57  J  -  l.MRK 

SUMSQ  -  SUMSQ  +  (DIST(J)  -  MEAN)**2 

57  CONTINUE 

STD  -  0.0 

IF(MRK  .GT.  2)  THEN 

STD  -  SQRT ( SUMSQ/REAL ( MRK  -  1)) 


52 


ENDIF 

PRINT*, 'MEAN  -  '.MEAN,  '  STD-  ',STD,'  MRK  -  ' ,MRK 
WRITE  (  7 1  ,*)  MEAN ,  STD ,  MRK 

CALL  SORT(DIST , ORDER_Y , START , MRK) 

CALL  SORT ( EVAR , ORDER_X , START , MRK) 

****  identifies  where  the  top  five  quality  messages  are  *** 
****  present  in  rank  order  of  the  estimated  quality  ******* 

DO  54  J  -  1,  COUNT 
DO  53  M  -  1,5 

IF(INT(ORDER_Y(J) ) . EQ. INT(ORDER_X(M) ))THEN 
IF(  BINO(I , J)  .NE.  1)THEN 
BINO(I , J)  -  1 
ENDIF 
ENDIF 

53  CONTINUE 

54  CONTINUE 


****  calculate  the  correlation  coefficient  *^***irk********* 
IF  (  MRK  .GT.  3) THEN 

CALL  S PEAR ( ORDER_X , ORDER_Y , CORR , T_STAT , START , COUNT ) 

PRINT*, STIME.CORR.T  STAT  ,MRK 

ENDIF 

WRITE(60,333)STIME,CORR,T  STAT, MRK 
333  FORMAT ( 3X , 14 , 5X , F8 . 5 , 2X , F10 . 5 , 2X , 15 ) 

STIME  -  STIME  +  DTIME 

DO  55  K  -  1, COUNT 
EVAR(K)  -  0.0 
DIST(K)  -  0.0 
ORDER_X(K)  -  REAL(K) 

ORDER_Y(K)  -  REAL(K) 

55  CONTINUE 

60  CONTINUE 

DO  70  I  -  1, SLICES 

WRITE(60, 61) (BINO(I , J)  ,J-1,25) 

61  F0RMAT(1X, 2512) 

70  CONTINUE 


STOP 

END 
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APPENDIX  3 


SUBROUTINE  FIX(COUNT, SLICES, DATA, EDATA , COVMTX) 

*  This  subroutine  opens  a  real  position  file,  and  finds  the 

*  position  in  the  input  file  where  the  simulator  did  not 

*  output  a  position.  It  then  writes  in  the  missing  position. 

*  For  the  estimated  position  and  covariance  matrix,  it  writes 

*  data  in  for  the  missing  data  from  the  pervious  simulation 

*  increment.  Since  only  a  few  points  failed  to  be  written, 

*  this  will  not  bias  the  results. 

*  by  Eric  R.  Bechhoefer,  SMC  1089,  NPS  MONTEREY,  CA  93943 

* 


IMPLICIT  NONE 


*******shared 


REAL* 8 
REAL*8 
REAL*8 

REAL*8 

REAL*8 

*******local 


DATA(3, 2500) !real  eci  positions,  output 
EDATA( 3, 2500) {estimated  eci  positions,  output 
COVMTX(3 , 3 , 2500) {estimated  eci  covariance  matrix, 

! output 

COUNT  {the  number  of  samples,  input 

SLICES  {the  number  of  time  increments, 

! input 

variable************************^**************** 


REAL*8  HOLD ( 3 )  {temporary  holding  variable 

REAL*8  MHOLD (3,3  )  '.temporary  holding  variable  for 

{covariance  matrix 

INTEGER  PLACE (100) {array  that  hold  the  position  of  error 
INTEGER  ERROR  {counter  for  error 

INTEGER  PTR  '.pointer 

INTEGER  MRK  ! counter 

INTEGER  I , J ,K,L,N,M  {counter 
************************************************************* 


*  Initialized  the  variables 


PTR  -  0 
MRK  -  1 
ERROR  -  0 
DO  15  I  -  1,  100 
PLACE ( I )  -  0 

15  CONTINUE 

DO  17  I  -  1,3 
HOLD ( I )  -  0.0 
DO  16  J  -  1,3 

MHOLD ( I , J )  -  0.0 

16  CONTINUE 

17  CONTINUE 


*  read  in  the  first  simulation  run  real  positions  *********** 

DO  20  I  -  1, SLICES 

READ(25 ,*)  DATA(1 , I) ,DATA(2,I) ,DATA(3,I) 

PTR  -  PTR  +  1 
20  CONTINUE 

*  test  these  positions  against  the  remaining  run  positions, 

*  note  the  position  error,  and  place  the  correct  point  in 

30  CONTINUE 

IF(  PTR  .LT.  SLICES*COUNT)  THEN 
PTR  -  PTR  +  1 

IF(MRK  .GT.  SLICES)  MRK  -  1 

READ (25,*) DATA ( 1 , PTR ) , DATA ( 2 , P TR ) , DATA ( 3 , PTR ) 
IF(DATA(1 , PTR)  .NE.  DATA(1 , MRK) ) THEN 
DO  40  J  -  1,3 

HOLD(J)  -  DATA(J.PTR) 

DATA(J.PTR)  -  DATA(J.MRK) 

DATA( J , PTR+1 )  -  HOLD(J) 

40  CONTINUE 

ERROR  -  ERROR  +  1 
PLACE (ERROR)  -  PTR 
MRK  -  MRK  +  1 
PTR  -  PTR  +  1 

PRINT*, 'ERROR  -  '.ERROR, 'PTR  -  '.PTR 
END  IF 

MRK  -  MRK  +  1 
GOTO  30 
ENDIF 


*  read  in  the  estimated  position  data  and  place  in  missing  data 

PTR  -  0 
MRK  -  1 
50  CONTINUE 

IF  (PTR  .LT.  SLICES*COUNT)  THEN 
PTR  -  PTR  +  1 

READ(20 ,*)EDATA(1 , PTR) , EDATA(2 , PTR) , EDATA(3 ,PTR) 

IF (MRK  .LE.  ERROR) THEN 

IF(PTR  .EQ.  PLACE(MRK) )  THEN 

PRINT*, 'FOUND  ERROR  HERE,  PTR  -  ’ ,PTR 
DO  60  J  -  1,3 

HOLD(J)  -  EDATA(J.PTR) 

EDATA(J.PTR)  -  EDATA(J , PTR  -  SLICES) 

EDATA(J , PTR  +  1)  -  HOLD(J) 

60  CONTINUE 

PTR  -  PTR  +  1 
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MRK  -  MRK  +  1 
ENDIF 
END  IF 
GOTO  50 
ENDIF 
PTR  -  0 
MRK  -  1 


*  read  In  and  fix  the  covariance  matrix* 

70  CONTINUE 

IF  (PTR  .LT.  SLICES*COUNT)  THEN 
PTR  -  PTR  +  1 

READ(50,*)COVMTX(1,1,PTR) ,C0VMTX(1,2,PTR) ,C0VMTX(1,3,PTR) 
READ(50 ,*)COVMTX(2 , 1, PTR) , COVMTX(2 , 2 , PTR) , COVMTX(2 , 3 , PTR) 
READ(50,*)C0VMTX(3 , 1 , PTR) ,COVMTX(3 , 2 , PTR) , COVMTX(3 , 3 , PTR) 

COVMTX( 2 , 1 , PTR)  -  COVMTX(l , 2 , PTR) 

COVMTX( 3 , 1 , PTR)  -  COVMTX(l, 3 ,PTR) 

COVMTX(3 , 2 , PTR)  -  C0VMTX(2 , 3 , PTR) 

IF (MRK  .LE.  ERROR) THEN 

IF(PTR  .EQ.  PLACE (MRK))  THEN 

PRINT* , ' FOUND  ERROR  HERE ,  PTR  -  ' , PTR 
DO  90  J  -  1,3 
DO  80  K  -  1,3 

MHOLD(J.K)  -  COVMTX(J ,K,PTR) 

COVMTX(J ,K, PTR)  -  COVMTX(J ,K,PTR  -  + 

SLICES) 

COVMTX(J ,K,PTR  +  1)  -  MHOLD(J.K) 

80  CONTINUE 

90  CONTINUE 

PTR  -  PTR  +  1 
MRK  -  MRK  +  1 
ENDIF 
ENDIF 
GOTO  70 
ENDIF 

RETURN 

END 
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APPENDIX  4 


r 


SUBROUTINE  FIX1 (COUNT , SLICES , DATA , EDATA , COVMTX) 
************************************************************** 

*  This  subroutine  is  similar  to  FIX,  however,  it  is  designed 

*  to  run  with  the  CROSS  program.  It  identifies  error  by  using 

*  a  error  file,  which  contains  the  real  eci  positions  at  every 

*  scan  period.  Values  that  can  be  identified  are  used  to  fill 

*  in  for  the  missing  data. 

*  by  Eric  R.  Bechhoefer,  SMC  1089,  NPS  Monterey,  Ca  93943 

* 

IMPLICIT  NONE 

******shared  variables**************************************** 
REAL  DATA(3,2500)  !real  eci  position 

REAL  EDATA(3 , 2500)  ! estimated  eci  position 

REAL  COVMTX(3 , 3 , 2500)  ! estimated  eci  covariance  matrix 

INTEGER  COUNT  ! the  number  of  sensors 

INTEGER  SLICES  !the  number  of  time  periods 


******local  variables***************************************** 


REAL  D(3 , 2500) 
REAL  E(3 , 2500) 

REAL  C(3 , 3 , 2500) 
REAL  X(3 ,48) 

INTEGER  PLACE (500) 


•raw  real  eci  position  data 
!raw  estimated  eci  position  data 
!raw  estimated  eci  covariance 
! array  of  the  real  eci  positions 
!  used  as  the  check 
! array  containing  error  position 
!  in  the  data  files 


INTEGER  ERROR 
INTEGER  PTR 
INTEGER  MRK.HK 
INTEGER  I , J ,K,L,N,M 


!how  many  errors  were  found 
!a  pointer 
! markers 
! counters 


******begin  code********************************************** 


OPEN (66 , FILE  -  ' ERROR ', STATUS  -  'old') 

PTR  -  1 
MRK  -  0 
ERROR  -  0 
DO  15  I  -  1,  100 
PLACE(I)  -  0 
15  CONTINUE 

**read  in  the  booster  position  for  flight  time**************** 
DO  20  I  -  1, SLICES 

READ(66 ,*)  x(l,I),x(2,I),x(3,I) 

20  CONTINUE 
MRK  -  1 

**read  in  the  raw,  eci  real  position  data  files**************** 
30  READ(25,*,END  -  33) (D(J , PTR) , J-l , 3) 
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PTR  -  PTR  +  1 
GOTO  30 
33  CONTINUE 

PTR  -  1 
HK  -  0 

40  IF(  hk  .LT.  SLICES*COUNT)  THEN 
HK  -  HK  +  1 

**check  for  deviation  of  data  file  position  from  key,  note** 

IFCMRK  .GT.  SLICES)  MRK  -  1 
IF(D(1 , PTR)  .EQ.  X(1 ,MRK) )  THEN 
DATA(l.HK)  -  D(l.PTR) 

DATA(2 ,HK)  -  D(2,PTR) 

DATA(3 ,HK)  -  D(3,PTR) 

PTR  -  PTR  +  1 

ELSE  !put  in  value  that  is  easy 

DATA(1,HK)  -  10.0  !to  identify  as  error 
DATA(2,HK)  -  10.0 
DATA(3 ,HK)  -  10.0 
ERROR  -  ERROR  +  1 
PLACE (ERROR)  -  HK 
ENDIF 

MRK  -  MRK+1 

GOTO  40 
ENDIF 

**finished  reading  in  eci  real  position  file,  identified 
**errors*************************************************** 
PTR  -  1 
MRK  -  1 

PRINT*, 'FINISHED  READING  IN  REAL  DATA' 

50  READ(20 ,* , END  -  55) (E(J , PTR) , J-l , 3) 

PTR  -  PTR  +  1 
GOTO  50 

**finished  reading  in  estimated  position  file*************** 
55  CONTINUE 

PTR  -  0 
HK  -  0 

PRINT*, 'READ  IN  ESTPOS','  ERROR  -  '.ERROR 
**reorder  data  to  include  missing  data********************* 
60  IF(PTR  .LT.  SLI CES*COUNT ) THEN 
PTR  -  PTR  +  1 
IF (MRK  .LE.  ERROR) THEN 

IF(PTR  .EQ.  PLACE(MRK) )  THEN 

PRINT*, 'FOUND  ERROR  HERE,  PTR  -  ' ,PTR 
DO  65  J  -  1,3 
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EDATA(J.PTR)  -  0.0 
65  CONTINUE 

MRK  -  MRX  +  1 
ELSE 

HK  -  HK  +  1 
EDATA(l.PTR)  -  E(1,HK) 

EDATA(2 , PTR)  -  E(2,HK) 

EDATA(3 , PTR)  -  E(3,HK) 

ENDIF 
ENDIF 
GOTO  60 
ENDIF 

**finished  reading  in  estimated  position,  start  reading**** 
**in  the  covariance  matrix  data**************************** 
PTR  -  1 
MRK  -  1 
HK  -  0 

70  READ (50,*, END  -  75) (C(l , J , PTR) , J  -  1,3) 

READ (50,*) (C(2,J ,PTR) ,J  -  1,3) 

READ(50,*,END  -  75) (C(3 , J , PTR) , J  -  1,3) 

PTR  -  PTR  +  1 
GOTO  70 
75  CONTINUE 

PRINT*, 'READ  IN  THE  COV  MTX' 

PTR  -  0 

80  IF  (PTR  .LT.  SLICES*COUNT)  THEN 
PTR  -  PTR  +  1 
IF (MRK  .LE.  ERROR) THEN 

IF(PTR  .EQ.  PLACE (MRK) )  THEN 
DO  90  J  -  1,3 
DO  85  K  -  1,3 

COVMTX(J,K,PTR)  -  1.0 
85  CONTINUE 

90  CONTINUE 

MRK  -  MRK  +  1 
ELSE 

HK  -  HK  +  1 
DO  110  J  -  1,3 
DO  100  K  -  1,3 

COVMTX(J,K,PTR)  -  C(J,K,HK) 

100  CONTINUE 

110  CONTINUE 

ENDIF 
ENDIF 
GOTO  80 
ENDIF 

RETURN 

END 
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APPENDIX  5 


SUBROUTINE  SORT  (X,  ORDER,  START,  COUNT) 
************************************************************** 

*  This  subroutine  takes  an  array,  bubble  sorts  it  in  ascending 

*  order  and  returns  an  array  ORDER  that  holds  that  order 

*  By  Eric  R.  Bechhoefer,  SMC  1089,  NPS  Monterey,  Ca  93943 

IMPLICIT  NONE 

*shared  variables********************************************* 


REAL* 8  X(100)  ! the 

! ordered,  input 

REAL* 8  ORDER(IOO)  ! the 

INTEGER  START  ! 

INTEGER  COUNT  ! the 


data  that  needs  to  be 

order  of  the  data,  output 
number  of  time  increments 


*local  variabl e s********************************************** 


REAL* 8  HOLD 
REAL*8  HOLD_A 
INTEGER  FIRST 

INTEGER  LAST 
INTEGER  J 
LOGICAL  SORTED 


!  temporary  holding 
.'temporary  holding  for  order 
! start  sorting  at  this  part  of 
!  the  array 

!  sort  the  array  to  this  point 
! if  sorted  then  true 


*  START  OF  CODE  ********************************************* 
SORTED  -  .FALSE. 

FIRST  -  START 
LAST  -  START  +  COUNT  -  2 
5  IF( .NOT.  SORTED)  THEN 
SORTED  -  .TRUE. 

DO  10  J  -FIRST,  LAST 

IF(X(J) .GT.X(J+1) )  THEN 
HOLD  -  X(J) 

HOLD_A  -  ORDER(J) 

X(J)  -  X(J+1) 

ORDER(J)  -  0RDER(J+1) 

X(J+1)  -  HOLD 
0RDER(J+1)  -  HOLD_A 
SORTED  -  .FALSE. 

ENDIF 

10  CONTINUE 

LAST  -  LAST  -  1 
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GOTO  5 
END  IF 
RETURN 
END 


APPENDIX  6 


SUBROUTINE  SPEAR  (X , Y , R_SUB_S , T_STAT , START, COUNT) 
******************************************************************** 

*  This  subroutine  takes  in  two  arrays  that  contain  rank  order, 

*  and  calculates  the  rank  correlation  coefficient  and 

*  t-statistic  associated  with  it.  It  assumes  that  there  are 

*  rela  ively  few  ties  in  the  rank  order. 

*  by  E'ic  R.  Bechhoefer,  SMC  1089  NPS ,  Monterey,  Ca  93943 

* 

IMPLICIT  NONE 

♦shared  variables******************************************** 

REAL*8  X(2500)  !an  array  containing  rank  order 

REAL*8  Y(2500)  !an  array  containing  rank  order 

REAL*8  R_SUB_S  ! calculated  rank  correlation,  out  REAL*8 

T_STAT  ! calculated  t-statistic,  out 

INTEGER  START  !ptr  indicating  the  start  of  the 

’array  for  calculations 

INTEGER  COUNT  ! the  number  of  time  increments 

♦local  var i ab 1 e s ********************************************** 

REAL* 8  DSQAR(3000)  .'array  that  hold  the  difference 

! squared  of  X  and  Y 

REAL*8  SUN_DSQAR  ! the  sum  of  DSQAR 

INTEGER  LAST  ! calculate  the  correlation 

(through  last 

INTEGER  I 

*  BEGIN  CODING  *********************************************** 

*  initialize  variables**************************************** 

LAST  -  START  +  COUNT  -  1 
SUM_DSQAR  -  0.0 
DO  5  I  -  START, LAST 
DSQAR(I)  -  0.0 
5  CONTINUE 

*  calculate  the  sum  of  squares******************************** 

DO  10  I  -  START,  LAST 

DSQAR(I)  -  (X(I)  -  Y(I) )**2 
10  CONTINUE 

DO  20  I  -  START,  LAST 
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SUM_DSQAR  -  SUM_DSQAR  +  DSQAR(I) 
20  CONTINUE 


*  calculate  the  correlation  coefficient  ********************* 

R_SUB_S  -  1.0  -  ( 6*SUM_DSQAR) / 

+  REAL (COUNT  *  (COUNT**2  -  1)) 

IF(R_SUB_S  .GT.  .995)  R_SUB_S  -  .995 

*  calculate  the  t— statistic  *********************************** 

T_STAT  -  R_SUB_S  *  SQRT(  REAL (COUNT  -  2))  / 

+  SQRT ( 1  -  R_SUB  S**2) 


RETURN 

END 


APPENDIX  7 


PROGRAM  MIX 

*  This  program  takes  an  ECI  position  of  a  sensor  platform  and 

*  generates  25  sensor  position  centered  on  this  point,  varied 

*  by  2  degrees 

*  by  Eric  R.  Bechhoefer,  SMC  1089  Monterey,  Ca  93943 

* 

IMPLICIT  NONE 

*  local  variables  ******************************************* 

REAL* 8  RH0(2)  !  array  that  holds  the  rho  for 

!  position  and  velocity 
REAL*8  PHE(2)  !  phe  for  velocity  and 

!  acceleration 
REAL*8  THETA(2) 

REAL*8  THETAD(2,10)  !  transformed  data 

REAL*8  PHED(2,10)  !  transformed  phe 

REAL* 8  X(2) ,Y(2) ,Z(2)  !  eci  xyz  for  position  and  accel. 

REAL*8  PI  !  constant 

PARAMETER  (PI  -  3.1459265359) 

REAL* 8  COORD(2 , 3 , 50)  !  transformed  eci  coordinates 

REAL* 8  MKR 

INTEGER  PTR , I , J , K , M 

OPEN (30 ,  FILE  -  '/  POS  DATA', STATUS  -  'NEW') 

*  initialize  variables  ****************************** 

DO  10  I  -  1,10 
DO  5  J  -  1,2 

PHED(J.I)  -  0.0 
THETAD(J.I)  -  0.0 
5  CONTINUE 

10  CONTINUE 

DO  30  I  -  1,50 
DO  20  J  -  1,3 
DO  15  K  -  1,2 

COORD (K,J , I)  -  0.0 
15  CONTINUE 

20  CONTINUE 

30  CONTINUE  * 

*  ENTER  THE  ECI  COORDINATES  X,Y,Z  for  position  and 

*  acceleration  *  1 

X(l)  - 
Y(l)  - 
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Z(1)  - 
X<2)  - 
Y(2)  - 
Z(2)  - 

HKR  -  2.0  *  PI  /  180.0 


*  start  the  transformation 


DO  40  M  -  1,2 

RHO(M)  -  SQRT(X(M)**2  +  Y(M)**2  +  Z(M)**2) 

PHE(M)  -  DACOS (  Z(M) /RHO(M)  ) 

THETA(M)  -  DASIN  (Y(M) /  (RHO(M)  *  DSIN(PHE(M) ) ) ) 


PHED(M, 1) 
PHED(M, 2) 
PHED(M, 3) 
PHED(M,4) 
PHED (M , 5 ) 


PHE(M)  -  2*MKR 
PHE(M)  -  MKR 
PHE(M) 

PHE(M)  +  MKR 
PHE(M)  +  2*MKR 


THETAD(M, 1)  -  THETA(M) 
THETAD(M, 2)  -  THETA(M) 
THETAD (M , 3 )  -  THETA (M) 
THETAD(M,4)  -  THETA(M) 
THETAD (M, 5)  -  THETA (M) 


-  2*MKR 

-  MKR 

+  MKR 
+  2*MKR 


40  CONTINUE 
PTR  -  1 


DO  70  I  -  1,5 
DO  60  J  -  1,5 
DO  55  M  -  1,2 


COORD (M,l, PTR)  -  RHO(M) 

*  DSIN(PHED(M,I)) 

+ 

*  DCOS( THETAD (M,J)) 

COORD (M, 2, PTR)  -  RHO(M) 

*  DSIN(PHED(M, I) ) 

+ 

*  DS IN( THETAD (M, I)) 

COORD (M, 3, PTR)  -  RHO(M) 

*  DCOS(PHED(M, I) ) 

PRINT*, (COORD(M,K, PTR) , 

K-  1,3) 

WRITE(30 ,*) (COORD(K, PTR) ,  K-  1,3) 

55 

CONTINUE 

PTR  -  PTR  +  1 

60 

CONTINUE 

70 

CONTINUE 

STOP 

END 
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